/*=============================================================================
*
*   Copyright (C) 2020 All rights reserved.
*
*   Filename:		basic_arithmetic.hpp
*
*   Author: Wang Zhecheng - wangzhecheng@yeah.net
*
*   Date: 2020-01-17 16:53
*
*   Last Editors: Wang Zhecheng - wangzhecheng@yeah.net
*
*   Last modified:	2020-07-24 18:07
*
*   Description:
*
=============================================================================*/

#ifndef BASIC_ARITHMETIC_H
#define BASIC_ARITHMETIC_H

#include "../types.hpp"
#include <limits>

namespace geoxfem {
namespace arithmetic {

template <typename T, typename P>
bool approximate(const T& a, const P& b, const T& tol = 1e-3) {
  auto diff = std::abs(a - b);
  if(diff > std::abs(a) * tol && diff > std::abs(b) * tol && diff > tol)
    return false;
  return true;
}

template <typename T>
inline int sgn(T val) {
  return (T(0) < val) - (val <= T(0));
}

inline _real sind(const _real& _degree) {
  return _real(
    std::sin(fmod((_degree), _real(360)) * _real(M_PI) / _real(180)));
}
inline _real cosd(const _real& _degree) {
  return _real(
    std::cos(fmod((_degree), _real(360)) * _real(M_PI) / _real(180)));
}
inline _real tand(const _real& _degree) {
  return _real(
    std::tan(fmod((_degree), _real(180)) * _real(M_PI) / _real(180)));
}

/*
 * Constexpr version of the square root
 * Return value:
 *	- For a finite and non-negative value of "x", returns an approximation
 *for the square root of "x"
 *   - Otherwise, returns NaN
 */
constexpr _real sqrt_iter(_real x, _real curr, _real prev) {
  return curr == prev ? curr : sqrt_iter(x, r_5 * (curr + x / curr), curr);
}
constexpr _real sqrt(_real x) {
  return x >= 0 && x < std::numeric_limits<_real>::infinity()
           ? sqrt_iter(x, x, 0)
           : std::numeric_limits<_real>::quiet_NaN();
}
inline int pow(const int& _a, const size_t& n) {
  int ans = 1;
  for(size_t i = 0; i != n; ++i)
    ans *= _a;
  return ans;
}
}  // namespace arithmetic

}  // namespace geoxfem

constexpr _real root2           = geoxfem::arithmetic::sqrt(2.);
constexpr _real root3           = geoxfem::arithmetic::sqrt(3.);
constexpr _real root5           = geoxfem::arithmetic::sqrt(5.);
constexpr _real root2Reciprocal = 1. / root2;
constexpr _real root3Reciprocal = 1. / root3;
constexpr _real root5Reciprocal = 1. / root5;

#endif /* BASIC_ARITHMETIC_H */
